Coulomb and Liquid Dimer Models in Three Dimensions 
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We study classical hard-core dimer models on three-dimensional lattices using analytical ap- 
proaches and Monte Carlo simulations. On the bipartite cubic lattice, a local gauge field generaliza- 
tion of the height representation used on the square lattice predicts that the dimers are in a critical 
Coulomb phase with algebraic, dipolar, correlations, in excellent agreement with our large-scale 
Monte Carlo simulations. The non-bipartite FCC and Fisher lattices lack such a representation, 
and we find that these models have both confined and exponentially deconfined but no critical 
phases. We conjecture that extended critical phases are realized only on bipartite lattices, even in 
higher dimensions. 



The statistical mechanics of dimers on a lattice that 
interact with one another only via hard-core exclusion 
has long been of interest to mathematicians and physi- 
cists [1, 2]. It is one of the simplest models describing 
the arrangement of anisotropic objects on a regular sub- 
strate. Applications [3] range from diatomic molecules 
on surfaces to spin ice in a magnetic field [4] . 

By Kasteleyn's theorem, on two-dimensional (2d) pla- 
nar lattices, the statistical mechanics of (close-packed) 
dimer coverings can be computed exactly [5]. A consis- 
tent picture has emerged from this work for a large class 
of 2d dimer models. On bipartite 2d lattices, dimer mod- 
els are in confined phases in which the free energy of two 
inserted test monomers (unpaired sites) increases with 
separation. The increase is logarithmic for phases with 
algebraic dimer correlations and linear for the remain- 
ing ones. An example of the former is the square lattice 
[5, 6] and one of the latter is the exotic diamond-octagon 
"4-8" lattice which exhibits two confining phases as the 
strength of the diamond bonds is varied [3]. By con- 
trast, non-bipartite lattices exhibit both deconfined and 
confined phases but always with exponentially decaying 
dimer correlations except at the boundary between such 
phases. Examples are the triangular [7, 8] and kagome 
[9] lattices, which are deconfined with exponentially de- 
caying correlations, and the Fisher lattice - equivalent to 
the 2d Ising model [10] - which exhibits a deconfinement 
transition [11]. 

Dimer models on two-dimensional bipartite lattices can 
also be understood through their height representations 
[12, 13]. Within this powerful framework, the two sub- 
categories of critical and non-critical dimer correlations 
are described as "rough" and "flat" phases. In cither case 
the defect interaction corresponding to the monomer free 
energy is long-ranged. Dimer models on non-bipartite 
lattices lack a similar long-wavelength description. 

In this paper we generalize the height representation 
to three dimensions. We show that dimer models on 
bipartite 3d lattices admit a local gauge representation 
which results in a "Coulomb" phase with algebraic, dipo- 



lar forms for the dimer correlations and monomer inter- 
actions that fall off inversely with their separation. Wc 
present large scale Monte Carlo simulations on 3d lat- 
tices and demonstrate that the dimer model on the cubic 
lattice is critical. On non-bipartite lattices we show that 
there exist both confined and deconfined phases, with 
exponentially decaying dimer correlations and monomer 
interactions, finding that the face-centered cubic (FCC) 
lattice and the 3d Fisher lattice realize these cases. As 
will be discussed elsewhere, the Coulomb and liquid 
phases can be identified with the U(l) and deconfined 
phases [14] in the corresponding quantum dimer models 
[15]- 

We first report our results on the cubic lattice along 
with the gauge representation. Thereafter we consider 
the FCC and Fisher lattices and conclude with a sum- 
mary and some conjectures on related models. 
Cubic lattice: Generally, for any dimer configuration, 
we define dimer numbers n(x) as ji,(x) = 1 if the bond 
between sites x and x+e^ is occupied, and zero otherwise 
(ej is the unit vector in direction i). Close-packed hard- 
core dimers obey the condition ^ Jij(x) — e,;) = 1. 
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FIG. 1: Dimer configuration on the bipartite square lattice 
(left) and corresponding configuration of the divergence- free 
magnetic field (right — a dark arrow carries 3 times the flux 
of a light one). The two lower quadrants represent plaque- 
ttes which can be flipped (dimers rotated by n/2). These 
plaquettes have zero average B-field. 

On a bipartite lattice, each dimer touches one site on 
each sublattice (cf Fig. 1, where the two sublattices are 
indicated by light and dark dots). Using a sublattice sign 
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factor cr x = ±1 depending on whether the site x belongs 
to one or the other sublattice, we can now define a field 
variable on each link 

Bi(x) ^CTxMx)-;?- 1 ) (1) 

(where z is the coordination number) as shown in Fig. 1. 
For close-packed dimers the lattice divergence of the field 
B vanishes, Xli-^i( x ) — B t (x — e,) = i.e., B is a lat- 
tice magnetic field without monopoles. A monomer is a 
monopole with charge ±1 depending on the sublattice. 

We define the lattice flux t^s through a surface E which 
does not contain any sites as the sum of the magnetic 
fields on the links piercing it 0s = (JJ B-dS)i at tice- For a 
cube with periodic boundary conditions, the flux through 
any surface that wraps around the system is invariant 
under local rearrangements of the dimers and under lat- 
tice translations of the surface. In particular, if we let 
Si be planes perpendicular to the cubic unit vectors e*, 
the fluxes <3?j through them are the maximal invariants 
that characterize a given topological sector of the dimcr 
model. For an L 3 cube the maximal possible flux is L 2 /2. 

A dimer configuration may be represented by a lattice 
magnetic field in any dimension. In 2d, one solves the 
constraint V • B = through B=Vx/i [16] where h, the 
height function mentioned earlier, is a scalar field on the 
dual lattice. In 3d, B is computed from a vector potential 
A on the links of the dual lattice. Given an arbitrary A 
we recover B as its (lattice) curl, B = V x A which is 
computed as 

W=E^ ( 2 ) 

yedp 

the oriented sum of the link variables on the boundary dp 
of the dual plaquette Pi(x) pierced by the link (x, x+e,). 
While in 2d, h is defined up to a global constant, in 3d, 
we have local gauge transformations, Aj(y) — > Aj(y) + 
A(y + ej) — A(y), where A is an arbitrary function on the 
sites of the dual lattice. Consequently, it is necessary to 
pick a gauge to work out properties of the A. The fluxes 
can be computed from the lattice line integrals of the 
vector potential along the boundaries of £j. Hence the 
sector with all $j = is obtained from gauge fields A 
that obey periodic boundary conditions themselves. 

In analogy to 2d we now conjecture that the long wave- 
length fluctuations of A and therefore B are governed by 
a probability distribution for the coarse-grained fields: 

-K [ (VxA) 2 -f f B 2 

P[A\ cx e 2 Jv y ' = e 2 Jv (3) 

in the vicinity of the zero flux <E>i = state. In the expo- 
nent of Eq. (3) the energy B 2 /2 of a magnetic field ap- 
pears naturally. Configurations that locally minimize the 
(coarse-grained) field strength (cf the lower quadrants of 
Fig. I) maximize the number of flippablc plaquettes with 
two parallel dimers and have high entropy, as described 
by Eq. 3. 



Two comments are in order: 

i) The assertion Eq. 3 implies that the gauge field is in a 
Coulomb phase, in the language of lattice gauge theories. 
The existence of this phase in our lattice system is not 
in conflict with Polyakov's proof of confinement [17] for 
the standard U(l) lattice gauge theory in 3d because in 
our case the microscopies explicitly forbid the monopoles 
that were crucial to his analysis. 

it) Gauge invariance explictly forbids any relevant op- 
erators at the fixed point defined by Eq. 3; this is the 
standard explanation of the masslessness of the photon. 
Consequently the prediction of a Coulomb phase is self- 
consistent and weak perturbations cannot give rise to 
anything new. This should be contrasted with the sit- 
uation in 2d, where vertex operators can become rele- 
vant even at weak coupling and, depending on the height 
stiffness and the radius of the height field, lead to a flat 
phase instead of the rough phase described by a purely 
Gaussian action. 

Returning to the ansatz (3) it is straightforward to 
deduce the long distance correlator, 

^) B M) = ^ XiXj ~ rHij (4) 

which is of the standard 3d dipole form. Reinserting the 
sublattice sign factors gives the connected dimer correla- 
tors. 

To test the dipole form, we have carried out Monte 
Carlo simulations using the pocket algorithm [7, 19] on 
large cubic lattices of size L 3 with L up to 128, with 
periodic boundary conditions. We have computed the 
connected correlation function for dimers at site x and at 
x+r, both pointing in the same direction ei = [100]. The 
vector r was taken as a multiple of lattice vectors [111] 
[010], [100], [110]. The correlations nicely fall off as l/r c 
in the regime o « r « i, with £ close to 3. Furthermore, 
the ratios of the correlations also agree very well with the 
predicted dipole form Eq. (4), as shown in Fig. 2. 

Going beyond this regime, we can even compute the 
scaling form of the correlations for r <~ L through a 
proper treatment of the periodic boundary conditions 
used by our algorithm, which explores all topological 
sectors. In the long wavelength description, we can 
write the fields in a given flux (topological) sector as 
B(x) = £\ ($ 4 /£ 2 )e 4 + B'(x) so that the field B' now 
carries zero flux. When $ = is nonzero, the 

stiffnesses K\\ and Kj_ for fluctuations parallel and per- 
pendicular to it will no longer be equal. Nevertheless 
these are semi-microscopic and hence can only depend 
on the average local magnetic field. Consequently we ex- 
pect that both K l± = K + 0($ 2 /L 4 ) where K is the 
stiffness in the zero flux sector. 

With these considerations we can generalize (3) to 

P[B] cx e 2 l e 2 Jv y 11 H ± ' . (5) 
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FIG. 2: Monte Carlo data of connected correlations (sym- 
bols) between parallel dimers, z 2 (ni(x)ni(O)) — 1, in various 
directions, multiplied by the sublattice sign factor a x and by 
R 3 , plotted vs. the Euclidean distance between the dimers, 
R. The coordination z — 6 for the cubic lattice. Left: For a 
system of size L = 128, 7? <C L, the dipolar asymptotes are in- 
dicated by horizontal lines. The overall scale factor is fixed by 
normalising those lines with respect to the [100] correlations. 
Right: For a system size L = 32, the data are compared to 
our finite-size formula; the only fitting parameter is an overall 
scale factor. 



The first factor will ensure that <& ~ 0(VZ) whence 
K^j_ = K + 0(1/ L 3 ) so that the stiffness anisotropy 
can be ignored for large systems. 

We can now deduce the correlations in a given flux 
sector, (B,(x)B J (0)) # = Q^j/L 4 + (fl<(x)B<(0)) where 
the second piece is independent of <I>. After averaging 
the first term with weight e~"?^, it equals Sij / (KL 3 ). 
From dimensional analysis we know that the second piece 
is of the form fij(x/L)/(KL 3 ) so that the two terms 
together have the appropriate form for a finite size scaling 
function. 

Finally, we have also determined the correlator 
(£^(x)£>j (0)) in a finite-sized sample subject to periodic 
boundary conditions. We are not aware of a closed form 
evaluation of this quantity but can write it as a sum in 
momentum space, which needs ultraviolet regularization. 
As we eventually compare it to a lattice simulation, it is 
most convenient to do this via a lattice sum The result- 
ing forms, e.g. (B x (ma)B x (0)) = C/(KL 3 ) with C given 

by 

[4 — 2 cos(27r ariy/L) — 2 cos(27r an z /L)] i2* an m 
+ £? [6-2£ iC os(27ran 4 /L)] & ' ' 



can then be compared directly with the simulations' re- 
sults for distances large compared to a lattice constant, 
a. In Fig. 2, the results are shown for a system of size 
L = 32. After adjusting the one free parameter to 
Ka 3 = 4.9 the curves agree well with simulation data 
for r larger than a few lattice spacings. A more detailed 
finite-size analysis that lets the exponent £ vary finds it 



to be C = 3.00 ± 0.02, in excellent agreement with the 
dipolar form. 

The gauge representation can be used to compute other 
operators of interest. Most notably, it shows that the 
interaction between two monomers is an attractive en- 
tropic force with the same form as that between oppo- 
sitely charged monopoles, i.e. an inverse squared force. 
One can also consider the Wilson loop (exp(i J c A ■ dl)) 
along a closed contour C. In the deconfined phase this 
will exhibit a perimeter law at large loop sizes, i.e. it 
will decay as e~ p where P is the perimeter of the loop 
C, although with interesting corrections coming from the 
long range nature of the Coluomb force [18]. 
FCC Lattice: As was emphasized above, we expect 
a fundamental distinction between bipartite and non- 
bipartite dimcr models. For the latter, the lack of a gauge 
representation indicates the absence of a Coulomb phase. 
To investigate this, we have simulated dimers on the face- 
centred cubic lattice, a simple non-bipartite 3d Bravais 
lattice. 

In Fig. 3, we display the connected correlations be- 
tween two parallel dimers in two different directions on a 
system of size L = 10, i.e. containing 4 x 10 3 sites. The 
decay of the oscillating correlations is extremely rapid, 
and fits well to an exponential form with a correlation 
length of £ = 0.35 ± 0.01 nearest neighbour distances. 
This establishes that the FCC dimcr model is in an ex- 
ponentially deconfined phase. 
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FIG. 3: Monte Carlo data for normalized connected dimer 
correlations, |144(rii(x)ni(0)> — 1| on the L = 10 FCC lattice, 
for i along the [110] direction. Separations x are along [110] 
and [011] and are measured in nearest-neighbour spacings, as 
indicated. 

Fisher Lattice: It is also possible to identify a 3c? dimer 
model with exponential correlations and both confined 
and deconfined phases analytically. This is done by map- 
ping the 3c? Ising magnet on the cubic lattice onto a dimcr 
model on the decorated cubic lattice shown in Fig. 4, 
where the original bonds have fugacity z = 1/ tanhif 
(K being the strength of the ferromagnetic coupling) 
whereas those internal to the decorating clusters have fu- 
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gacity 1. The correlator of the cubic Ising model between 
spins Si, Sj at sites i and j can be expressed in terms of 
monomer correlators of the resulting dimer model in the 
same way as described in detail in Ref. 11 for the cor- 
responding 2c? lattice [10] - in fact, the analysis carries 
over to any dimension. 




FIG. 4: The 3d Fisher lattice is obtained by decorating the 
sites of a simple cubic lattice (left) with the cluster shown on 
the right. Its six external legs (heavy lines) correspond to the 
original bonds of the cubic lattice. The numbers denote the 
possible locations of the test monomers. 

One obtains (SiSj) = Sfe^fc^i m (i,fe»),0',fej) ; where 
m (i,ki)-U-kj) denotes the correlator of test monomers lo- 
cated on sites ki of cluster i and kj of cluster j (see Fig. 4) . 
This correlator is given by the ratio of the partition func- 
tions with and without the pair of monomers present. 

In the high temperature phase of the Ising model, 
(SiSj) — > exponentially, so that all monomer pairs 
are confined. At low temperatures (SiSj) decays ex- 
ponentially to a constant. Therefore at least one pair 
of monomers is deconfined. Unless possible algebraic 
terms present in the four correlators exactly cancel, which 
seems very unlikely, the deconfinement is exponential. 
Conclusions: We have argued that the distinction be- 
tween the behaviour of bipartite and non-bipartite lat- 
tices familiar from dimer models in 2c! holds also in 3d. 
The former are characterized by a conservation law in a 
long wavelength description which gives us a Coulomb 
phase in 3d and mutatis mutandis should do so in all 
d > 3 as well. Non-bipartite lattices lack such a con- 
straint and can be generically expected to exhibit expo- 
nentially confined or deconfined phases as they do in our 
3d examples. The Coulomb phase should also arise in 
other problems that involve conservation laws, e.g. the 
statistical mechanics of spin ice which is also the Ising 
antifcrromagnet on the pyrochlore lattice [20]. Finally, 
we note that Hcrmclc, Fisher and Balents [21] are con- 
sidering a set of closely related three-dimensional models. 
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